##Copy the part starting with source() into the console and run it there - this will create a log file that includes the output

#source("~/work/November_2024/R_Script_Refinement_Figures.R", echo = T, max.deparse.length=10000)
rm(list=ls(all=TRUE))
sink("~/work/November_2024/R_out/log_files/full_log_refinementPLOTS.log", split = T)

# ======================================================================
# Project:    Closed borders, closed minds? COVID-related border closures,
#             EU support and hostility towards immigrants
#
# Script:     Appendix Figure E1 
#
# Authors:    Lisa Herbig (l.j.herbig@uva.nl)
#             Asli Unan (a.unan@uva.nl)
#
# Date:       2025-03-24
# ======================================================================

# Description:
# This script generates Appendix Figure E1, presenting covariate balancing analysis for the border closures model.

# ----------------------------------------------------------------------
# Load libraries
# ----------------------------------------------------------------------

# Covariate Balancing Plots 
library(reshape2)
## Read the files in
df_iso_work_dd = readRDS("~/work/November_2024/R_out/df_iso_work_dd_25.rds")

# EU

EU_N <- PanelMatch(lag = 4, 
                   time.id = "weekyear_i", 
                   unit.id = "iso2.y_i", 
                   matching=FALSE, 
                   treatment = "any_closure_binary_control0", 
                   refinement.method = "none",  
                   covs.formula = ~ East +working_abroad_n+ kr_foreigner_n +CONTACTABROAD_n+ HOMEOFFICE_n +mode_phone_n+ kr_inz_rate_n+age_n + female_n + unemployed_n + university_n + EU_residence_n + iso2_pop_n + Voting_2018_LeftRight_n, 
                   data = as.data.frame(df_iso_work_dd), 
                   match.missing = FALSE,
                   size.match = 1,
                   qoi = "att", 
                   lead = 0:4, 
                   outcome.var = "con_europe_n", 
                   placebo.test = T)

cvb_EU_1 <- as.data.frame(get_covariate_balance(EU_N$att,as.data.frame(df_iso_work_dd), covariates = c("East", "working_abroad_n" ,"kr_foreigner_n", "CONTACTABROAD_n" ,"HOMEOFFICE_n" ,"mode_phone_n" , "kr_inz_rate_n", "age_n", "female_n" , "unemployed_n" , "university_n" , "EU_residence_n" , "iso2_pop_n", "Voting_2018_LeftRight_n")))


EU_m <- PanelMatch(lag = 4, 
                   time.id = "weekyear_i", 
                   unit.id = "iso2.y_i", 
                   matching=TRUE, 
                   treatment = "any_closure_binary_control0", 
                   refinement.method = "mahalanobis",  
                   covs.formula = ~ East +working_abroad_n+ kr_foreigner_n +CONTACTABROAD_n+ HOMEOFFICE_n +mode_phone_n+ kr_inz_rate_n+age_n + female_n + unemployed_n + university_n + EU_residence_n + iso2_pop_n+ Voting_2018_LeftRight_n, 
                   data = as.data.frame(df_iso_work_dd),
                   match.missing = TRUE,
                   size.match = 1,
                   qoi = "att", 
                   lead = 0:4, 
                   outcome.var = "con_europe_n", 
                   placebo.test = T)

cvb_EU_2 <- as.data.frame(get_covariate_balance(EU_m$att,as.data.frame(df_iso_work_dd), covariates = c("East", "working_abroad_n" ,"kr_foreigner_n", "CONTACTABROAD_n" ,"HOMEOFFICE_n" ,"mode_phone_n" , "kr_inz_rate_n", "age_n", "female_n" , "unemployed_n" , "university_n" , "EU_residence_n" ,"iso2_pop_n", "Voting_2018_LeftRight_n")))

cvb_EU_1 <-apply(as.data.frame(apply(cvb_EU_1, 1, abs)),2, mean)
cvb_EU_2 <-apply(as.data.frame(apply(cvb_EU_2, 1, abs)),2, mean)

cvb_plot <- as.data.frame(rbind(cvb_EU_1, cvb_EU_2))

cvb_plot$Refinement <- factor(c("None","Mahalanobis"), levels= c("None", "Mahalanobis"))

cvb_plot <- melt(cvb_plot, id.vars=c("Refinement"))
colnames(cvb_plot)[2] <- "Time Period"
colnames(cvb_plot)[3] <- "Average of Mean Differences"

ggplot(cvb_plot, aes(`Time Period`, `Average of Mean Differences`,  color = Refinement, group = Refinement)) +
  geom_hline(yintercept=0, color = "black", size=0.9)+
  geom_point(size=4.5, shape=18) +
  geom_line(size=1)+
  theme_bw(base_size = 22)

ggsave(path = "~/work/November_2024/R_out/refinement_plots/", filename = "EU_RefinementPlot.pdf")



# Refugees

IMI_N <- PanelMatch(lag = 4, 
                    time.id = "weekyear_i", 
                    unit.id = "iso2.y_i", 
                    matching=FALSE, 
                    treatment = "any_closure_binary_control0", 
                    refinement.method = "none",  
                    covs.formula = ~ East +working_abroad_n+ kr_foreigner_n +CONTACTABROAD_n+ HOMEOFFICE_n +mode_phone_n+ kr_inz_rate_n+age_n + female_n + unemployed_n + university_n + EU_residence_n + iso2_pop_n + Voting_2018_LeftRight_n, 
                    data = as.data.frame(df_iso_work_dd), 
                    match.missing = FALSE,
                    size.match = 1,
                    qoi = "att", 
                    lead = 0:4, 
                    outcome.var = "immigration_pca_n", 
                    placebo.test = T)

cvb_IMI_1 <- as.data.frame(get_covariate_balance(IMI_N$att,as.data.frame(df_iso_work_dd), covariates = c("East", "working_abroad_n" ,"kr_foreigner_n", "CONTACTABROAD_n" ,"HOMEOFFICE_n" ,"mode_phone_n" , "kr_inz_rate_n", "age_n", "female_n" , "unemployed_n" , "university_n" , "EU_residence_n" , "iso2_pop_n", "Voting_2018_LeftRight_n")))


IMI_m <- PanelMatch(lag = 4, 
                    time.id = "weekyear_i", 
                    unit.id = "iso2.y_i", 
                    matching=TRUE, 
                    treatment = "any_closure_binary_control0", 
                    refinement.method = "mahalanobis",  
                    covs.formula = ~ East +working_abroad_n+ kr_foreigner_n +CONTACTABROAD_n+ HOMEOFFICE_n +mode_phone_n+ kr_inz_rate_n+age_n + female_n + unemployed_n + university_n + EU_residence_n + iso2_pop_n+ Voting_2018_LeftRight_n, 
                    data = as.data.frame(df_iso_work_dd),
                    match.missing = TRUE,
                    size.match = 1,
                    qoi = "att", 
                    lead = 0:4, 
                    outcome.var = "immigration_pca_n", 
                    placebo.test = T)

cvb_IMI_2 <- as.data.frame(get_covariate_balance(IMI_m$att,as.data.frame(df_iso_work_dd), covariates = c("East", "working_abroad_n" ,"kr_foreigner_n", "CONTACTABROAD_n" ,"HOMEOFFICE_n" ,"mode_phone_n" , "kr_inz_rate_n", "age_n", "female_n" , "unemployed_n" , "university_n" , "EU_residence_n" ,"iso2_pop_n", "Voting_2018_LeftRight_n")))

cvb_IMI_1 <-apply(as.data.frame(apply(cvb_IMI_1, 1, abs)),2, mean)
cvb_IMI_2 <-apply(as.data.frame(apply(cvb_IMI_2, 1, abs)),2, mean)



cvb_plot <- as.data.frame(rbind(cvb_IMI_1, cvb_IMI_2))

cvb_plot$Refinement <- factor(c("None","Mahalanobis"), levels= c("None", "Mahalanobis"))

cvb_plot <- melt(cvb_plot, id.vars=c("Refinement"))
colnames(cvb_plot)[2] <- "Time Period"
colnames(cvb_plot)[3] <- "Average of Mean Differences"

ggplot(cvb_plot, aes(`Time Period`, `Average of Mean Differences`,  color = Refinement, group = Refinement)) +
  geom_hline(yintercept=0, color = "black", size=0.9)+
  geom_point(size=4.5, shape=18) +
  geom_line(size=1)+
  theme_bw(base_size = 22)


ggsave(path = "~/work/November_2024/R_out/refinement_plots/", filename = "Refugee_RefinementPlot.pdf")

